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Abstract 

A point mass at the center of an ellipsoidal homogeneous fluid is 
used as a simple model to study the effect of rotation on the shape and 
external gravitational field of planets and stars. Maclaurin's analytical 
result for a homogenous body is generalized to this model. The absence 
of a third order term in the Taylor expansion of the Maclaurin function 
leads to further simple but very accurate analytical results connecting 
the three observables: oblateness (e), gravitational quadrupole (J2), 
and angular velocity parameter (q). These are compared to observa- 
tional data for the planets. The moments of inertia of the planets are 
calculated and compared to published values. The oblateness of the 
Sun is estimated. Oscillations near equilibrium are studied within the 
model. 
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1 Introduction 



The rotation induced oblateness of astronomical bodies is a classical problem 
in Newtonian and celestial mechanics (for the early history, see Todhunter 
^). It has twice played an important role in the history of science. In 
the early eighteenth century measurements indicated a prolate shape of the 
Earth, in strong conflict with the Newtonian prediction. This was later 
shown to be wrong by more careful measurements by Maupertuis, Clairaut, 
and Celsius in northern Sweden in 1736. Then, in 1967, measurements of the 
solar oblateness were published and according to these it was much larger 
than the Sun's surface angular velocity would explain. The confirmation 
of general relativity by Mercury's perihelion precession would then be lost. 
Also this problem is now gone and the modern consensus is that the solar 
oblateness is too small to affect this classic test of general relativity jH El 15 • 

The subject of the flattening of rotating astronomical bodies is thus 
quite mature. The classical theory is due mainly to Clairaut, Laplace, and 
Lyapunov. Also Radau, Darwin, de Sitter and many others have made im- 
portant contributions. More recent accounts of the theory can be found in, 
for example, Jeffreys Zharkov et al. 6_, Cook ^7 , Moritz j8j and, partly, 
in Chandrasekhar [0]. Some pedagogical efforts can be found in Murray 
and Dermott JHI) or in Kaula II. As is plain from these references the 
theory is quite involved. Only the unrealistic assumption that the body is 
homogeneous gives compact analytical results. Otherwise a specified radial 
density distribution is needed and one resorts to cumbersome series (mul- 
tipole) expansions, or purely numerical methods, for quantitative results. 
Here we will present analytical results based on the assumption that the 
body consists of a central point mass surrounded by a homogeneous fluid, 
the so called point core model. By varying the relative mass of the fluid and 
the central point particle one can interpolate between the extreme limits of 
Newton's homogeneous body and the Roche model with a dominating 
small heavy center. The point core model goes back to the work of G. H. 
Darwin. More recently it has been used to study the shape of outer planet 
moons, see Hubbard and Anderson ^31) Dermott and Thomas jl4j . 

Apart from the basic point core approximation (i) several further approx- 
imations are assumed here. These are: (ii) that the shape is determined by 
hydrostatic equilibrium and (iii) that the shape is ellipsoidal. These are not 
consistent. According to Hamy's theorem, see Moritz |H], the exact shape 
is not ellipsoidal, so we regard an ellipsoidal shape as a constraint and flnd 
the equilibrium shape, among these, by minimizing the energy. A further 
approximation (iv) is the neglect of differential rotation. On the other hand 
fixed volume ("bulk incompressibility" ) need not be assumed; the equilib- 
rium volume problem separates from the shape problem. In spite of these 
approximations, which are standard in the literature, the mathematics can 
be quite involved. In this article I hope to clarify and simplify it as much 
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as possible. 

Mathematically our model then becomes a three degree of freedom me- 
chanical system for which we can calculate the kinetic and potential (grav- 
itational plus centrifugal) energies exactly. Multipole expansions in terms 
of spherical harmonics are not needed. The three degrees of freedom cor- 
respond to the three semi-axes of the ellipsoid (a, b, c) but these are trans- 
formed to three generalized coordinates that describe size (or volume), R, 
spheroidal, ^, and triaxial, r, shape changes, see Eq. ©. The statics problem 
of equilibrium shape is solved by minimizing the potential energy, U{R, ^, r) 
for a fixed R, Eq. (|47j) . For slow rotation the shape will not be triaxial so 
r = 0. Finding the shape, or flattening, is then only a matter of minimizing 
a dimensionless potential energy. 



given by the last two terms of Eq. (|48|) . Here is defined by Eq. (|45|) and 
plotted in Fig. ^ while is a constant that depends on rotational parameter 
q and dimensionless moment of inertia k^, Eq. (|55|) . The root £,{k) of the 
equation V''(C) — = thus gives the flattening. This root is to very high 
accuracy given by 



see Eqs. (jl)T|) and (|T^ . for the rotational parameters that can be found in 
the solar system. 

Prom this result, simple formulas, Eqs. (|(i2|) - (|64jl . relating the observ- 
ables, rotation parameter, q, gravitational quadrupole, J2, and excentricity 
squared, e^, are obtained. These appear to be, partly, new, and their useful- 
ness is demonstrated by comparing with empirical data for the Sun and the 
rotating planets of the solar system. Variational methods have been used 
before to study similar problems, see for example Abad et al. [T5]. Denis 
et al. have pointed out that variational methods generally fail to pro- 
vide estimates of their accuracy. Therefore the agreement of our formulas 
with empirical data, as demonstrated in Table /cite(table), is important 
and demonstrates that our model catches the essential physics of rotational 
flattening. 

Finally small amplitude oscillations near the equilibrium are investi- 
gated, starting from the Lagrangian 



for our three degree of freedom model system. This gives useful insight 
into the physics of free stellar or planetary oscillations and their coupling to 
rotation. Our approximations are, however, too severe for these results to 
be of quantitative interest. 



(1) 





L = T{R,^,T,R,tf)-U{R,tT) 



(3) 
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2 Basic geometric quantities 



Assume that x, y, z are rectangular Cartesian coordinates in three-dimen- 
sional space of a point with position vector 



r = xe^ + ysy + ze^. 
An ellipsoid, with semi-axes a, b, c, is the solid defined by 



-I +if) +1- 



1 1/2 



< 1. 



If we put 
and define 



a = (e + t)R, b=iC- t)R, c = (e^ - t 



2\-l 



R 



+ 



y 



+ 



the ellipsoid is the set of points that fulfill 

E{r;tT)<R, 

where R is the geometric mean (or volumetric) radius R 
formula for the volume of an ellipsoid now gives 

V = ATTahc/2, = inR^/S, 



(4) 

(5) 
(6) 

(7) 
(8) 



^abc. The 



1/2 



(9) 



and we see that changes of ^ and r do not change the volume, only the shape 
of the ellipsoid. 

In what follows we will also consider the special case of spheroids. A 
spheroid is an ellipsoid with two equal semi-axes. We take a = b ^ c. This 
means that we take r = in the formulas above so that 



^R, c = r'R, 



and (O becomes 



+ 



,-| 1/2 



A spheroid is thus defined as the set of points r that fulfill 

E{r;(„0) < R. 



(10) 



:ii) 



(12) 



Here R is the (geometric) mean radius and ^ is a parameter that determines 
the shape in such a way that ^ < 1 corresponds to a prolate spheroid, ^ = 1 
to a sphere, and ^ > 1 to an oblate, or flattened, spheroid. 
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More familiar parameters used to define the shape of a spheroid are the 
eUipticity e and the eccentricity e. The eUipticity is defined by 

e^^ = i-r^. (13) 

This is sometimes also called the (geometric) oblateness or the flattening 
(denoted /). Solving this equation, and expanding around e = 0, we have 

One notes that e is positive for oblate and negative for prolate spheroids 
respectively. The eccentricity, e > 0, is defined as the usual eccentricity of 
the ellipse that is the intersection of the spheroid and a plane containing the 
z-axis. One finds that 

e^-V^ = i-r^ (15) 

gives the eccentricity of these ellipses in the oblate case. For the prolate 
case = 1 — 



3 The mass distribution 

Consider a non-rotating spherically symmetric body consisting of point par- 
ticles with masses rrii and position vectors 

rO = x^e^ + + ^°e,. (16) 

Assume that rotation induces a deviation of the positions so that the new 
positions, = f{rf), are given by 

ri{^, t) = Xie:c + Viey + Zie^ = {i + r)xie^ + (^ - T)y°ej; -I- (^^ - T^)"^2;°e^. 

(17) 

This implies that we assume that the elastic displacement field can be 
parameterized by the shape parameters ^, r with non-rotating positions 
corresponding to ^ = l,r = 0. One notes that points r° that obeyed 
E{r^; 1,0) = r, i.e. were found on a sphere of radius r < R, move to the 
surface given by 

E{r;^,r)=r. (18) 

The sphere is thus assumed to deform to an ellipsoid. 

Assume that a body initially is spherically symmetric and has a mass 
density goir), where r is the ordinary distance from the center. It is tempting 
to assume that a natural deformation of the body when it starts to rotate is 
to a spheroidal shape in such a way that the originally spherical equidensity 
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surfaces, go{r) = constant, deform to similar spheroidal surfaces given by 
-E(r;^,0) = r. Unfortunately this is an approximation for real bodies. For 
a density that increases towards the center, the ellipticity of the equidensity 
surfaces decrease with r, in a manner described by Clairaut's equation, see 
[HI Ei IZl Note, however, that for a constant density, or a constant 

density with a central point mass, the assumption is not approximate, since 
any surface is an equidensity surface, only the shape of the outer surface 
matters. 

Thus we now assume such a mass distribution. To be precise we assume 
that a mass Mi is distributed with constant density 

^1 = Mi/y (19) 

inside the surface, E{r] ^, r) = R, and that the point mass Mq is at the 
center, r = 0, of the ellipsoid. The density g is thus given by 

1 U otherwise 
where 5-? is the three-dimensional Dirac delta function. 



4 Moments of inertia and quadrupole tensor 

If we put rrii for the masses of the particles, the inertia tensor of the body is 
given by 1^(3 = J2i mix^iX-yi Sap-XaiXfn), where we have put, Xi = xu, yi = 
X2i, Zi = X3i, and 6af3 is the usual Kronecker delta. The inertia tensor for 
the non-rotating body is diagonal with all diagonal elements equal and given 

by 

/ = Mi^Rf = II = Y^m,[{x^f + (21) 

i 

Here we have introduced the total mass, 

M = Y,mi = Mo + Mi, (22) 
i 

as well as the dimensionless radius of gyration k. The square of this pa- 
rameter, K?, also appears in the literature denoted I /MR^ (or C/Ma?', or 
simply k) and is also called the moment of inertia coefficient. 

We calculate the dimensionless moment of inertia k? for the undeformed 
Q of Eq. ((201). Using Eq. ((STJ and r'^ = x^ + y"^ + we find 

f 2 3iVfi 2 

/ = MURf = / Qiix"^ + y'^)dV = - / —^r'^AiTr'^dr = -MiR^, 
Jr<R 3 Jo 4:7: R-^ 5 

(23) 
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where dV = Airr'^dr and Eq. (|19j) have been used. This is the usual moment 
of inertia for a sohd sphere of mass Mi . This imphes that 

5 M 5 M ^ ' 

is the squared (dimensionless) radius of gyration of the undeformed body. 
Note that < < 2/5 since the model can go between the limits of a point 
(k^ — > 0) and a homogeneous sphere (k^ = 2/5 = 0.4). 
For the deformed body one gets, using (|17|) . 

Izz = J2M^l + yf) = MiKRfif + T^) = lie + r^). (25) 

i 

Izz is often denoted by C. For spheroids (r = 0) the other moments of 
inertia become Ixx = lyy = ^-^(1 + C~^)^^- For these one frequently finds the 
notation A, and respectively, in the literature. The quadrupole tensor, 
D^p , is given in terms of the inertia tensor by D^p = I^^ Sa/s — "^lap ■ It is 
identically zero for the undeformed body and diagonal for the same axes as 
the inertia tensor. For the spheroidal (r = 0) body its components are given 

by 



1 



D^D,., = Dyy = -I{l-rie, (26) 

and Dzz = —2D. Thus D = Izz — Ixx (or D = C — A) is positive for oblate 
and negative for prolate shapes. 



5 Potential and hydrostatic equilibrium 



One measure of the oblate shape that is of great interest is the, so called, 
gravitational quadrupole J2. This dimensionless quantity, gives the devi- 
ation of the gravitational potential, 0, outside the body, from the simple 
1/r-dependence of spherically symmetric bodies. The first two terms in the 
multipole expansion of this potential are, see e.g. Stacey HI 



GM 



J2 - 



a\2 Scos^T?- 1 



(27) 



where 'd is the polar angle (colatitude) in spherical coordinates, M the total 
mass, and a the equatorial radius of the body. One can show ^Hl that J2 
of formula is given by J2 = D/{Mo?). Use of D given in and 
a = ^R, then gives 



J2 



D 
Ma? 



1 2 2 



(28) 



for the quadrupole. Here we have also used a = (,R and Eq. ((T3)) . Note that 
J2 = <^2(e^) is independent of R, but directly proportional to the eccentricity 
squared. 
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Assuming the oblateness to be small, and due to hydrostatic equilibrium 
in the combined gravitational and centrifugal force fields, one can derive a 
formula connecting J2 with the angular velocity and the ellipticity. One uses 
that the shape of the surface is determined by the hydrostatic equilibrium 
equation 

Vp + qV{(P + (Pc) = 0, (29) 
where p is pressure, g mass density, and 

(f)^(r, T?) = sin^ d (30) 

is the potential of the centrifugal force in a system rotating with angular 
velocity to about the z-axis. This means that the surface of the body must 
be an equipotential surface of (/> + (/>c. The constancy of (/)(r, -d) + (jidf, '&) on 
the surface of the body can be used to derive the approximate relationship 

^=1^2 + \q (31) 

[Jj- Here we have introduced the ratio q of equatorial centrifugal to gravi- 
tational acceleration, 

^ " ~ GM/R^ ~ 'GM ~ 4^" ^ ' 

Here g = M/V is the mean density. Instead of q the notation m occurs in 
the literature, see Zharkov et al. 0, who reserve q for the corresponding 
quantity with the equatorial radius a replacing the mean radius R (recall 
R? = a^c). The advantage of the definition (|.S2|) is that it makes q indepen- 
dent of flattening. The relationship pif) . is only valid to first order in small 
quantities, so here the difference does not matter, but we will go to higher 
order below, so take notice. This relationship between the three quantities 
q (or m), e (or /), and J2, can thus be used to test the hypothesis of hydro- 
static equilibrium empirically. This will be done below, using more exact 
results. For Jupiter and Saturn these will improve significantly on (|^T|l . 



6 Rotational and gravitational energy 

To study a rotating body, we go to the rotating system, and assume that all 
particles are affected by the centrifugal and gravitational potential energies. 
If the system is rotating rigidly about the z-axis with angular velocity vector 
uj = uiez, the potential energy of particle i in the centrifugal force field is 
mi(f)c{ri,i}i). The total centrifugal potential energy of the body (system of 
particles) is 

^c = -^E"^^(^i +y'V' = -^^..w'- (33) 

i 
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This is the negative of the rotational kinetic energy. Using (|25j) we can write 

= + r^)^' = -\M{nRf{e + r^)^' (34) 

for the total potential energy of the particles in the centrifugal potential. 

The gravitational potential from an ellipsoid ©, of constant density qi, 
at an interior point is |19j 

/•oo / x'^ \ ds 

(t>{x, y, z) = -Ggiirabc / 1 ^ ^— (35) 



where P{s) = \/ (a^ + s){b'^ + s){c'^ + s). Since MQ(j){0, 0, 0) is the interac- 
tion energy of the homogeneous ellipsoid with a point particle of mass Mq 
at the origin one finds the expression 

^o--GM„M,-l ^ (36) 

for this interaction energy. Here abc = and Eq. (|19|) has been used. 

An elementary calculation based on see Landau and Lifshitz PO] , 
also shows that the gravitational self energy of the ellipsoid is 



Use of the definitions © and the substitution s = R^u gives 
ds 1 /""^ dti 



/o P(s) i?7o v/[(C + + u\ [(e - ry + u\ ae - r^)-^ + nj 
If we define 



(38) 



1 f°° ds 



we find that the total gravitational energy of the mass distribution (|20|1 is 
given by 

^ ^ ^ GMi /3Mi 3Mo\ , 

$ = $1 + $0 = (^—i + —0 j r), (40) 

where x(l)O) = —1. We define the constant ly, the dimensionless gravita- 
tional radius, through 

^ = -^x(C,r). (41) 

Comparing with (|4()j) and using M = Mq + Mi together with the result (|24j) , 
a small calculation reveals that the dimensionless gravitational radius, u, is 
given by 

l/i/ = (15/8)^2(2 -3^2). (42) 
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It is interesting to note that the gravitational energy of a density (|2fl|) . with 
fixed mass M and radius R, is minimized for = 1/3 (this is incidentahy 
very close Earth's value). 

Using the definition (|39|) it is easy to find the Taylor expansion of x 
around ^ = 1, r = 0. This gives 

x(e,r) = -l + ^(C-l)' + ^r2 + ... (43) 

so the sphere really corresponds to a minimum in the gravitational energy. 
Higher order terms are easily generated by computer algebra programs, but 
we will not need them here. To study the behavior further away from the 
minimum one might use, instead of ()39() . an excellent Pade approximation 
of the ellipsoidal energy derived by Dankova and Rosensteel |21j . Near the 
spherical minimum it is almost indistinguishable from the exact expression. 
For spheroids we define 

V'(0=X(C,0) (44) 

to be the dimensionless gravitational potential energy. For r = the integral 
in ()39|) can be found in terms of elementary functions. Several different 
alternative expressions exist and they are always different depending on 
whether ^ is greater or less than unity. However, noting that, arctanh \z = 
iarctanz, the two expressions in 



,^^arctanh\/l — 



arctan 



for < ^ < 1 (prolate case) 

(45) 

for 1 < ^ < cxD (oblate case). 



are equivalent and can be used for all < ^ < oo. These are really the same 
real function if only the correct branch is chosen when passing through ^ = 1. 
Series expansion, of either alternative, around ^ = 1 gives 

4 Q2 2Q 

m = -1 + 5 - 1)' - - 1)' + - 1)' + - l)% (46) 

This function, which is shown in Figure Q has the following properties. 
The prolate limit, ^ 0, of an infinite line is zero: limg^oV'(C) = 0- The 
oblate limit, ^ oo, an infinite circular disc, is also zero: lim^^oo V'('?) = 0- 
The sphere (^ = 1) corresponds to a minimum of V'(C) so that ^(1) = — 1, 
V''(l) = 0, and ^"(1) = 8/5. 

7 Minimizing the energy 

The total energy, U, of a static body (star or planet) will, apart from the 
centrifugal energy, $c of Eq. and gravitational energy, $ of Eq. ||1T|). 
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also consist of some volume dependent energy, ^y{R). We can then write 
this total energy it in the form 

C 1 

U{R,^,t) = <^,{R) + -^Xitr) - -M{^R)W{e + T^). (47) 

By minimizing this energy we can find out whether a system tends to be 
triaxial or spheroidal. For small angular velocities the spheroidal case is 
known to be relevant [HII^- We thus put r = here and get 

C M'^ 1 

U{R, = MR) + - -M{^R)We. (48) 

Using this we can find out if the body will be prolate or oblate and how this 
is affected by rotation. Since we concentrate on stars and planets surface 
forces and other neglected shape dependent forces are normally very small 
compared to the gravitational and centrifugal energies. 
The equilibrium equations can be written 

^ = ^^'(^) - MiKRfu'^ = 0, (50) 

where we have denoted differentiation by with a prime. If we multiply 
the first equation by R and the second by ^ and subtract we get, after 
rearrangement , 

d$„ GM"^ r , /AN lU^^^i 

-^ = i^m)+^'m- (51) 

This is essentially an equation of hydrostatic equilibrium. Note that it de- 
termines the equilibrium i?-value with no direct reference to the angular 
velocity w. The second equation, (fSO]) . can be written 

-^i^'iO = M{kR)W^, (52) 

and this expresses a balance of the shape dependent part of the gravitational 
force with the centrifugal force. If we put 

^ ^ ~GM^'^^ " ^^^^ 
this turns into the concise expression 

V''(e) = K- (54) 

Since ■0'(1) = this will cause a shift of the shape from the spherical equi- 
librium, at ^ = 1, to ^ > 1, as long as a; > 0. Use of shows that 

/c = (8/15) g/(2 - 3^2) (55) 

The equilibrium ^ is thus entirely determined by q and the moment of inertia 

^2 
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8 Finding the position of the minimum 



To investigate Eq. (|54|) we use the second variant for "0(0 ^^1- (ESI) and 
get, after some algebra, 



I) 



(3 + e - 1) arctan( _ i) _ 3^^6 _ i 



(56) 



This looks simpler if we introduce the new variable ?7(^) = \/^^ — 1. Use of 
Eq. (fT5|) shows that 



T] 



- 1 = e/Vl-e^ 



(57) 



so that r] can also be regarded as a function of the eccentricity, r]{e). Now 
we get 



k 



(3 + 7? ) arctan 7] — St] /rj 



(58) 



For the homogeneous case k = 2g/3 and this formula, written in terms of 
the eccentricity by means of the identity, arctan(e/Vl — e^) = arcsine, is 
equivalent to the celebrated Maclaurin's formula '9 from 1742. 

Taylor expansion of the Maclaurin function around e = only contains 
positive even powers of e. The coefficient of vanishes so a very good 
approximation for this function is implied by 



(3 + ry^) arctanr/ — 37/ /rf 



15 



+ 0(e« 



(59) 



The two term approximation is nearly perfect for all planetary e-values, see 
Figure 121 Saturn has maximum e = 0.43 and the absolute error at this value 
is 2 • 10~^ and the relative error is 4 • 10~^. The relative error even remains 
below 1% up to e = 0.75. Using this approximation we must solve 



15 ^ 7 



and solving this for gives 



7 
2 



l + —k 



(60) 



(61) 



Use of Eq. (|28() gives k? = 2J2/e^. Insertion of this into H55|) eliminates k? 
from k = vK^q. Some algebra then leads to the equation 

q{e^J2) = {l + e^/7)[e^ -3J2) (62) 

for q. This equation can be solved for J2 or for to give 

1 



J2{q,e^) 
e\q,J2) 




(63) 
(64) 
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respectively. These results, together with 



e(2 -e)^ 2e, 



(65) 
(66) 



e 



1- Vl-e2, 



are thus more exact versions of Eq. Already the approximation 46^/15, 
to the Maclaurin function, gives only a 3% error for Saturn, and simpler 
formulae, but use of these more exact versions removes all numerical con- 
siderations from the problem. 

9 Theory versus observational data 

Table ^ compares modern data with the predictions of Eqs. ()62|) - H66|) . 
Observed data for q, e, and J2 from Lodders and Fegley ^] are given in the 
first (horizontal) row of data for each of the major rotating planets. In the 
same line calculated from = e(2 — e) is then given. The last value in the 
row is the dimensionless moment of inertia = I /ME? as given by Lodders 
and Fegley [52]. These values are based on calculation. For original sources 
the reader may consult The Planetary Scientist's Companion by Lodders 



For each pair of the values q, J2, and e^, in the first row, the third value 
is calculated using - H64() and the results are displayed in the second 
row of data for each planet in Tabled e is then found from equation (|66)) . 
For Jupiter and Saturn there is a third row of data in Table Q This row 
is analogous to the second row except that the data are calculated from 
observational pairs of data using the first order formula l|31j) so that q is 
given by observational e and J2 according to q = 2e — 3J2, e is found directly 
from (|3H) using the observational q and J2 and so on. These rows illustrate 
the fact that the non-linear formulas of this work are more accurate than 
first order results, in spite of the point core approximation. Even for the 
Earth the first order formula gives discrepancies in the third digit. 

Tabled shows clearly that the observational data and theory are in rea- 
sonable agreement for the planets with the exception of Mars and Uranus. 
The perfect agreement for the Earth is probably due to the fact that the 
oblateness is determined from the shape of the geoid (the equipotential sur- 
face at sea- level). For Mars the disagreement may depend on the fact that 
oblateness is determined from geometrical shape rather than from the shape 
of an equipotential surface. Also, of course, the assumption of hydrostatic 
equilibrium may not be perfect for Mars. The problems with the Uranus 
data, and to some extent those for Neptune, are probably due to observa- 
tional uncertainty. 



and Fegley 
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10 Oblateness and moment of inertia 



Instead of eliminating the moment of inertia (or mass ratio) of Eq. ()24() 
from the equations as done above one can retain it and try to use observed 
data to get information about . Series expansion of Eq. ()58|) with r] 
expressed in terms of appropriate variables gives 

" ^l-A^r\-^ ^0{i) (67) 



2 - 14 V2 - 3k 

for the flattening, or oblateness, expressed in terms of q and the dimen- 
sionless moment of inertia. For a homogeneous body = 2/5 so e = 
|(? + + 0((7^) and the old Newtonian result, e = |g, is obtained, to first 
order. In the opposite (Roche) limit of dominating central mass (k^ = 0) 
and one finds e = + ^q"^ + 0{q^). 

In a similar way we get from J2 of Eq. (|28() that 

J2 = ^^q+-(^^\ q^ + 0{q'') (68) 
2 - 3^2^ 21 V2 - 3^2; ' ^ ^ 

is the dimensionless quadrupole moment of the rotating body. For the ho- 
mogeneous (k;2 = 2/5) body this implies J2 = \q + sfs^'^ + 0{q^). In the 
opposite limit of ^ one finds that J2 goes to zero. 
Our results (|H7I) and (|H5|) . to first order, are 

^('^') = ^zW^' ^^^^ 

>/2(k2) = Y^^^ (70) 

and imply that J2 = en^. These may be compared to similar expressions 
derived by G. H. Darwin using the Radau equation 

^^'^'^ = 16 + 25(2^-3.2)2'^' (^1) 

J2(.^) = ^i^^^i^^,. (72) 
^ ' 48 + 75(2 -3.2)2^ ^ ' 

These two have the same values at k? = 2/5, the homogeneous sphere, 
and the same derivatives with respect to .2 at this point, as (|69|) and (|7n|) . 
respectively. The simple expressions (|^^ and ((7n|) . derived here, however, 
also have natural limits for = 0, in contrast to e.g. (|72|) . which becomes 
negative for .2 < 2/15. 

We now use Eq. ^ i.e. = 2J2/e^ and Eq. (jMl) for e^{q,J2) to 
calculate k? in terms of the observed q and J2. The last entry in the second 
row for each planet of Table ^ gives , 

.2(g,J2) = 2J2/e2(g,J2), (73) 
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as calculated from observational, first row, q and J2. Since q and J2 are 
easier to measure accurately than flattening this should give more reliable 
values than using the observed ellipticity directly. In this way one obtains 
information about the radial mass distribution in the interior of the planet 
from purely external data. The resulting K^-values are given last in the 
second line for each planet in Table ^ The data obtained in this way are 
compared with published and tabulated values for the planets (Lodders 
and Fegley |22J). For Earth and Mars the agreement is reassuring. For 
the outer planets there is naturally a lot of uncertainty but at least the 
increasing trend from a minimum at Saturn is a common feature. 

One problem with the outer planets, and with the Sun in particular, is 
to decide what angular velocity, u, should be used to find q. Since these 
bodies have differential rotation some average must be used. There is a 
well defined theoretical average angular velocity (Essen |2S1) but it is not 
directly observable. Once an average angular velocity has been selected for 
the Sun one can use the present theory and a theoretical K^-value (0.059) 
to estimate the solar oblateness and J2. The definition of average angular 
velocity in indicates that an angular velocity near the radius of gyration 
kR ~ 0.24i? should be appropriate. The further assumption that angular 
velocity is constant on coaxial cylinder surfaces (Kippenhahn and Weigert 
jl2j . Ulrich and Hawkins |21]) indicates that an angular velocity near the 
poles of the Sun rather than equatorial angular velocity is relevant. 

A handbook value [22j for the polar angular velocity of the Sun is = 
2.1 • 10-6 rad s'^. Use of this gives the g-value (1.15 • IQ-^) in Tableffl The 
precise numbers are not important here. What is important is that, because 
of the small k^, an angular velocity near the pole rather than an equatorial 
angular velocity is used. Otherwise the flattening will be exaggerated. Use 
of g = 1.15 • IQ-^ and k? = 0.059 from Lodders and Fegley (22j, together 
with formulas H69|) and H70|) . gives the results in the top row of Tabled that 
is, e ~ 6.3 • 10"^, and J2 ~ 3.7 • 10^''. These values are, at least, of the same 
order of magnitude as the currently best motivated values |31 1^ HU UHl 
which are roughly e ~ 9 • 10"^ and J2 ~ 2 • 10^^. 

11 Small oscillations near equilibrium 

We now wish to study the small amplitude motions of a gravitating rotating 
nearly spherical body. We wish to know how the radius R is affected by 
the rotation so we take R to denote the constant non-rotating equilibrium 
radius and let RX be the variable radius. A is thus a dimensionless variable 
which is equal to unity at the undeformed geometry. We parameterize the 
positions of the particles of the system as follows 

r,(A, r) = i?A[(e + r)x*e,. + (^ - r)y*e, + {f - r^h^e,]. (74) 
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The xl,y*,z* are then the non-rotating equihbrium dimensionless position 
coordinates of the particles. Equilibrium corresponds toA = .^ = l,T = 0, 
see Section |2j A is a dilatation (or radial pulsation) degree of freedom, 
while C 7^ 1 corresponds to spheroidal deformations and r 7^ to triaxial 
deformations. 

The velocity of the particles are 

^,(A, i f, A, r) = nm + r)x*e, + (^ - T)y*e, + (f - T^)-^z*e,\ (75) 
+Rm + r)xle, + f)y*ey - 2{ii - Tf){e - T^)-^z*e,]. (76) 

An elementary computation shows that the kinetic energy, T = ^ J2i "mivf , 
becomes 



r(A,e,f,A,e,r) = ^ 



2Ar(l + l)Af + A^(l + §W^ + A^ f 1 + ^ " ^A^^r 



Here 5 = - r^, and / = 2R^ J2i mxf = MR'^k^. 

The potential energy is given in Eq. (|T7j) . To get an explicit potential 
energy we must find some expression for the volume dependent (pressure) 
energy ^^{RX). As a model for this part of the energy we take the simple 
expression 

1 GM^ 1 

where we must take ?i > if the corresponding force is to balance gravity 
and prevent collapse. Putting these together we thus find the total potential 
energy 



(79) 

If we use the definition of k in Eq. (|53|) we find the expression 

for the potential energy of the system. 

Combining T of Eq. l|77|) with U to the Lagrangian L(A, ^, f. A, ^, r) = 
T—U we can get the dynamics of this three degree of freedom system exactly 
by finding and solving the Euler-Lagrange equations of the system. Here we 
will instead assume small oscillations and make a normal mode analysis. We 
first taylor expand U to quadratic terms around A = ^ = l,r = 0, and then 
solve for the zeros of the first derivatives of this quadratic to get the position 
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of the minimum. The result is 



^0 



1 + 

1 + 
0. 



k{8 + 5k) 



8(n -k)- 5{n + 3k)k 

5k{n + k) 
8{n-k) -5{n + 3k)k 



n 



We now introduce new variables, p, a through 

/) = A - Ao, CF = i-io]. 



(81) 

(82) 
(83) 

(84) 



Then, of course, p = \,a = ^. In the kinetic energy we replace A, ^, r, in the 
coefficients with the equilibrium positions Aq, Co; "^o and expand to first order 
in k. In the potential energy we make the same replacement in the quadratic 
taylor expansion, expand to first order in k and throw away constant terms. 
The result is that 



T 

and that 
U = 



MR^K^ 



+ + f2 + ^-^{pa - a') + ^(Sa^ + f'] 
2 In 



(85) 
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Only the degrees of freedom p and a are coupled and only when /c > 0. For 
A; = the spheroidal a-mode and the triaxial r-mode are degenerate (have 
the same frequency). The value n = 4/5 is special since for that value all 
three modes have the same frequency. 

For the p-mode pressure is restoring force in the contracting phase of the 
motion and gravity in the expanding phase, for the other two modes, o", r, 
gravity alone acts to restore the spherical minimum. If we put 



0^ 

"0 



8GM 
15i?3 



GM 



(87) 



Eq. ^ shows that /c = and Eq. (jSOI), with 

For n > 4/5 we then get the approximate eigen frequencies 



2e, that k = ^e. 



5n 



-ijj 



1^2 + ( 1 



1 + 



-1 

15ny 
16 
15n 



(89) 
(90) 
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to first order in k. From tliis one easily finds the following first order results 



T 



Tr 




(91) 



(92) 



(93) 



for the corresponding periods. 

The free radial oscillation mode for the Earth is known to have 
Tp = 20.45 min = 1228 s. This means that one can calculate n and find it 
to be n = 13.19 for the Earth. Since the equilibrium radius of the rotating 
Earth is RXq Eq. (|81|) . and a small calculation, now shows that Earth's 
mean radius is 863 m larger due to rotation compared to the non-rotating 
case. The fairly large n also shows that the g mode periods are essentially 
independent of n. For the Earth one finds that the r-mode has period 
Tr = 83.22 min. Due to the rotational splitting, which is given by 



the (T-mode is 17 seconds shorter. The longest observed free oscillation 
period of the Earth has period 53.8 min (see Udfas [23 )■ The main explana- 
tion for the discrepancy is probably that the model neglects elastic forces. 
Though not quantitatively reliable from this point of view the model has 
the advantages of showing in a simple way how rotational splitting of the 
modes arise and the order of magnitude of such splittings. 

12 Conclusions 

Some apparently new results relating to the classic theory of the figure of 
rotating bodies have been presented. The basic model, a point mass at 
the center of a homogeneous fluid, is characterized by their mass ratio, and 
interpolates between the limits of an ellipsoidal homogeneous fluid and a 
body dominated by a small central mass concentration. It allows simple 
analytic treatment but is still flexible enough to correctly describe the es- 
sential hydrostatics of real rotating planets as well as stars. Such models are 
always useful, especially when one wishes to analyze, compare, and under- 
stand large numbers of observational data. In spite of the simplicity there 
is no perturbation order to which the results are valid; the nonlinearity of 
the basic equations can be retained. To be more precise Eq. (|59l) shows 
that the formulas are valid to seventh order in the eccentricity (within the 
basic model with its simplified mass distribution). As demonstrated by the 
numerical experiments on Jupiter and Saturn data, this is essential. In fact 
Table n indicates that the geometric oblateness of the surface equipotential 



AT = Tr-Ta = 27re/ilo 



(94) 
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surface of Mars and Uranus determined from observed q and J2 values using 
(jnU probably are more reliable than current observational data. Finally the 
dynamics of the model reveals the essentials of the coupling and rotational 
splitting of the most basic free oscillations of a planet without recourse to 
expansion in spherical harmonics. 
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body 




e 


7 


e 


K 


Sun* 


1.15-10-5 


6.3-10-'' 


3.7-10-^ 


1.3-10-5 


0.059 


Earth^ 
Earth^ 


3.45-10-^ 
3.45-10-^ 


3.35-10-3 
3.35-10-3 


1.08-10-3 
1.08-10-3 


6.69-10-3 
6.69-10-3 


0.3307 
0.323 


Marsl 
Mars^ 


4.57-10-3 
7.04-10-3 


6.48-10-3 
5.24-10-3 


1.96-10-3 
2.78-10-3 


1.29-10-2 
1.04-10-2 


0.366 
0.375 


Jupiter 

Jupiter^ 

Jupiter^ 


8.34-10-^ 
8.29-10-2 
8.56-10-2 


6.49-10-2 
6.51-10-2 
6.37-10-2 


1.47-10-2 
1.45-10-2 
1.54-10-2 


1.26-10-1 

1.26- 10-1 

1.27- 10-1 


0.254 
0.233 


Saturn^ 
Saturn^ 
Saturn"^ 


1.40-10-1 
1.41 - 10— "'^ 
1.47-10-1 


9.80-10-2 

9.73-lU 

9.44-10-2 


1.63-10-2 
1.67-10—2 
1.87-10-2 


1.86-10-1 
1.85-10-1 
1.89-10-1 


0.210 
0.176 


Uranus^ 
Uranus^ 


2.89-10-2 
3.50-10-2 


2.29-10-2 
1.98-10-2 


3.52-10-3 
5.55-10-3 


4.53-10-2 
3.92-10-2 


0.225 
0.179 


Neptune^ 
Neptune'^ 


2.56-10-2 
2.34-10-2 


1.71-10-2 
1.82-10-2 


3.54-10-3 
2.80-10-3 


3.39-10-2 
3.61-10-2 


0.24 
0.196 



Table 1: Values of q = R^lo'^ /GM, of oblateness, e, and gravitational 
quadrupole, J2, ellipticity, e, squared and dimensionless moment of inertia 
= I/MR^. 

1. The first row for eacli planet gives literature P2| data. These are observational except 
^2 which are based on theoretical calculations. 

2. The second row for each planet gives the corresponding calculated values as given by 
formulas (1621 - 16611 in such a way that for each pair of observational q, J2 and e'^ the 
missing third is calculated. The moment of inertia is calculated from observational q and 
J2 using formula 12811 . = U'ijfP' , with e2(q, J2) from formula 1641 . 

3. The third row for Jupiter and Saturn gives data calculated in a similar way to that in 
the second row except that the first order formula I31I I has been used to find the third 
value from two observational values. 

4. For the Sun e, J2 and e'' have been calculated from an observation based g-value 
discussed in the text and a theoretical |22| . using formulas and = 2e 
respectively. 
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0^ 0.2 0^4 0.6 0^8 ^1 

e 

Figure 2: The Maclaurin function /(e) = [(3 + r7^)arctanr/ — 3r/]/r/^ with 
rj = e/Vl — compared to the approximations 4e^/15 (lower curve) and 
(4e^/15)(l + e^/7) (upper curve). 



23 



